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Using Ecological Niche Modeling for Biodiversity Conservation Guidance in the Western Podillya 
(Ukraine): Reptiles. Tytar, V., Sobolenko, L., Nekrasova, О. Mezhzherin, S. — Maximum entropy 
niche modeling was employed as a tool to assess potential habitat suitability for 10 reptile species and 
to map their potential distribution in the Western Podillya (Ukraine). We used climate, topography and 
human impact (assessed by the Human Footprint) as predictor variables. “Isothermality”, “temperature 
seasonality” and the “mean temperature of coldest month” were three most important factors in predicting 
habitat suitability and distribution. A profound contribution to the modeling has been displayed by the 
Human Footprint, meaning that human infrastructure may benefit reptile species. Areas have been 
distinguished that in the first place should be of interest to biodiversity conservationists targeting reptiles 
and maps summarizing predicted habitat suitability and species richness were produced for guiding 
conservation efforts. 
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Моделирование экологической ниши как инструмент для планирования мероприятий, 
направленных на сохранение биоразнообразия Западного Подолья (Украина): рептилии. 
Титар B., Соболенко JI., Некрасова O., Межжерин С. — Моделирование экологической ниши 
методом максимальной энтропии было использовано для оценки условий пребывания 10 видов 
рептилий и их распространения на территории Западного Подолья (Украина). Среди предикторов 
были использованы показатели климата, рельефа и антропогенного воздействия (оценивается 
по интегрированному индексу «человеческий след»). Среди важнейших факторов, которые 
определяют пригодность и распределение мест обитания рептилий, были «изотермичность», 
«сезонность температурного режима» и «средняя температура самого холодного месяца». 
Также существенный вклад вносит индекс «человеческого следа». Это может означать, что 
инфраструктура, созданная человеком, может формировать благоприятные условия для рептилий. 
Отмечены территории, которые в первую очередь должны представлять интерес для охраны 
рептилий, и созданы карты, которые обобщают прогнозируемую пригодность среды обитания и 
видовое богатство исследованного региона, что может быть использовано для целенаправленных 
природоохранных мероприятий. 


Ключевые слова: Maxent, моделированиеэкологической ниши, моделированиераспространения 
видов, рептилии, Украина. 


Introduction 


Even though population declines have been documented, the reptile fauna is not often taken into account 
as a conservation objective. Fragmentation and habitat loss threatens reptiles with extinction, and habitat 
specialists, such as endemic or rare species, are more vulnerable to changes in habitat gradients. Estimates 
indicate that 15-36 % of the world’s species of reptiles are threatened (Bohm et al., 2013). This is higher than 
for bird species (12 %) and around the same level with the percentage for all mammals (25 %) (Brandon et al., 
2005). 
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Protected areas have the potential to aid in the conservation of reptiles by protecting their remaining habi- 
tat and eliminating the threats that are causing their decline. The need to recognize the distribution of a species 
in a specific area should be considered baseline information for developing studies and management plans for 
its conservation. Unfortunately, insufficient basic ecological information and the actual distribution of their 
populations has become a problem for planning conservation strategies. 

For most of species, the possibility of occurrence in a certain area can be predicted by species distribution mod- 
els (SDMs). SDMs are becoming an important method and have been widely used (Hanspach et al., 2010; Johnson, 
Gillingham, 2008). SDMs are techniques that using the relationship between species occurrence and environmental 
conditions to model the geographical ranges of suitable-habitat for the certain species (Varela et al., 2011). 

However, with the increase of human activities and the development of human society, many anthropogenic im- 
pacts have been generated, which means the occurrence of species might not only relate to environmental conditions 
but also might rely on anthropogenic factors. Similarly as SDMs can connect species occurrence with environmental 
conditions, SDMs could also become a tool to investigate empirical relationships between species occurrence and 
anthropogenic factors when anthropogenic factors are involved in building the models. 

A variety of distribution modeling methods are available for predicting the potential geographical range 
of a species, however the performances of most species distribution modeling methods are poor when sample 
size is small (for instance, < 10). Under these circumstances Maxent (Phillips et al., 2006) stands out because it 
has been found to perform best among many different modeling methods (Elith et al., 2006) and may remain 
effective despite small sample sizes. Maxent is a maximum entropy based machine learning program that esti- 
mates the probability distribution for a species’ occurrence based on environmental constraints (Phillips et al., 
2006). It requires only species presence data (not absence) and environmental variable (continuous or categori- 
cal) layers for the study area. 

Significantly, this study is a local proposal that addresses the distribution of reptile species in the Western 
Podillya in Ukraine, so that more effective strategies can be accomplished. In particular, our objectives were to 
(1) identify the factors associated with (species) habitat distribution; (2) predict potential distributions of the 
species using known presence observations; and ultimately (3) produce a regional map of reptile species rich- 
ness for conservation guidance. 


Methods! 


Occurrence Data Collection and Processing 

For generating the occurrence data set used in the modeling we digitized published presence survey data 
from L. Sobolenko (2010), I. Dotsenko (2003); T. Hrynchyshyn (2008); T. Kotenko, O. Kukushkin, O. Zinenko 
(2008); I. Skiľ’s’kyi, N. Smirnov, Г. Khlus, L. Meleshchuk (2008); №. Smirnov, I. Skil’s’kyi (2008); Vikyrchak, 
N. Smirnov (2014) and unpublished reports from O. Vikyrchak, M. Drebet, and O. Nekrasova. Georeferencing 
(in OziExplorer v. 3.95.4 m) was accomplished for 297 point data obtained for 10 species: Anguis fragilis Lin- 
naeus, 1758 (22), Coronella austriaca Laurenti, 1768 (29), Emys orbicularis (Linnaeus, 1758) (15), Lacerta agilis 
Linnaeus, 1758 (45), Lacerta viridis (Laurenti, 1768) (70), Zootoca vivipara (Jacquin, 1787) (16), Natrix natrix 
(Linnaeus, 1758) (39), Natrix tessellata (Laurenti, 1768) (26), Vipera berus (Linnaeus, 1758) (6) and Zamenis 
longissimus (Laurenti, 1768) (29). These species make extensive use of both wetland and upland habitat and 
have been negatively impacted by habitat loss and degradation resulting from wetland drainage, agriculture, 
and urbanization. As a consequence, over much of their range, these species persist in small populations iso- 
lated from conspecifics by habitat fragmentation. 

Sources for taxonomy are from the IUCN Red List of Threatened Species, version 2014.3 (http://www. 
iucnredlist.org). 

The logistic probabilities provide a relative indication of the likelihood of occurrence by the species, 
but they do not define predicted occurrence in the binary, presence/absence manner typically required by 
managers. To establish which areas are climatically suitable for the species, we used the 10 percentile training 
presence logistic threshold. This has the effect of conservatively identifying a region of highest fit that does not 
allow outlying points of presence to expand the predicted area of occupancy beyond a core region. Binary maps 
of species occurrence were created by reclassifying the integer raster datasets, where all values above the given 
threshold were assigned a value of 1 and all other values were assigned a value of 0. 

Final versions of habitat suitability maps were considered to benefit by a simple smoothing filter as 
one of the means for coping with observation bias (Home et al., 2007) and generalizing raster outputs: 3 x 3 
neighborhood filtering implemented in DIVA GIS (www.diva-gis.org/) was applied for this purpose. 

Creating contour data layers from grid data layers was accomplished with SAGA 2.1.4_Win32, a system 
for automated geo-scientific analysis (http://www.saga-gis.org/). 


Environmental Data Collection and Processing 

Temperature, precipitation and solar radiation have been proved that have significant correlations with 
the distributions of reptile species and activities of reptile individuals (Nicholson et al., 2005). Topographical 
factors have been stated as necessary factors for SDMs (Guisan et al., 1999). The topography predictor variables 


! The study region, model building and evaluation have been described in a previous paper (Tytar et al., 2015). 
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selected were: topographical wetness index (TWI) and aspect. Aspect is a proxy for the amount of solar radiation 
on the ground surface, the TWI combines a measure of the upslope area and slope to predict the hydrology of 
a given location (Sorenson et al., 2005). Small values represent upper catenary positions (dry), and high values 
represent lower catenary (wet) positions. 

Bioclim climate data are derived from the monthly temperature and rainfall values (Hijmans et al., 2005). 
It generates more biologically meaningful variables about temperature and precipitation. Bioclim climate data 
contains 19 factors represent annual trends, seasonality and extreme or limiting environmental conditions. 

The Human Footprint (HF) data set (Sanderson et al., 2002) was taken as a measure of anthropogenic 
impact. The HF layer combines layers representing human population pressure (population density), human 
land use and infrastructure (built-up areas, night-time lights, land use/land cover), and human access (roads, 
railroads, navigable rivers). The resultant layer is normalized by biome to reflect the continuum of human 
influence on the natural environment. The HF score ranges from zero to 100. 

All environmental data layers were spatial resolution rasters (~1 km) with the same extent and cell 
alignment, arranged by the QDSM plugin in QGIS 2.6.1 freeware (www.qgis.org). 

Due to the high levels of correlations between many environmental variables, we filtered the initial 
variable set of 23 predictors based on the results of multi-colinearity analysis. by calculating the Variance 
Inflation Factor (VIF): VIF = 1/(1-R’), where В? is the coefficient of determination. VIF values higher than 10 
are considered to lead to problematic levels of multi-collinearity (Craney, Surles, 2002). Factors with VIF values 
that exceeded 10 were not used in the model. The 9 selected factors are shown in table 1. 


Table 1. Predictors, remaining after the multi-colinearity analysis 











Predictor VIF 
bio 3: Isothermality (bio 2/bio 7*) 2.0 
bio 4: Temperature Seasonality (STD** x 100) 285 
bio 6: Min Temperature of Coldest Month, °C 2.4 
bio 13: Precipitation of Wettest Month, mm 3.3 
bio 14: Precipitation of Driest Month, mm 4.1 
bio 19: Precipitation of Coldest Quarter, mm 5.3 
Aspect 1.0 
Topographical wetness index (TWI) 1.2 
Human Footprint (HF) 1.0 





“Bio 2: Mean Monthly Temperature Range; bio 7: Temperature Annual Range (°C*10). 
7 STD — standard deviation. 


Results and discussion 


Factors associated with (species) habitat distribution 

The modeling identified a number of environmental variables that had high contribu- 
tion to generating the potential distribution prediction of reptiles in the study area. As it 
is the case with reptiles in general, their large-scale distribution is mainly environmentally 
dependent, due to physiological features, their behavior and ultimately fitness is influenced 
directly by environmental conditions (Huey, 1982; Huey, Kingsolver, 1989). 

Based on percent contribution, of all 9 variables initially used in the modeling process, 
the most important ones and those that best explained the environmental requirements of 
the species, were isothermality (bio 3), temperature seasonality (bio 4), minimum tempera- 
ture of the coldest month (bio 6) and human impact, measured by the Human Footprint 
(HF); precipitation of the wettest month (bio 13) and precipitation of the coldest quarter 
(bio 19) were of importance only for Vipera berus (table 2). The percent contribution of 
precipitation of the driest month (bio 14), aspect and topographical wetness index (TWI) 
in all cases was below 10 %, so they were not included to table 2. 

Isothermality (bio 3) was a major contributor to models for all the 10 reptile species, 
particularly Zamenis longissimus (67.7 %) and Lacerta viridis (60.8 %). All considered in the 
study species displayed a similar bell-shaped response of increased suitability as isothermal- 
ity increases above a certain value, indicated by the response curves generated using only 
the corresponding variable. Response curves also gave an indication of the range under 
which the variable reaches its optimum suitability. For Zamenis longissimus, for instance, 
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Table 2. The percent contribution of environmental variables (factors) in predicting the species 
geographic distribution models 





























, Factors associated with (species) habitat distribution > 10 % 
Species | | | 
bio 3 bio 4 bio 6 bio 13 bio 19 HF 
Anguis fragilis 24.4 - 17.3 - - 36.4 
Coronella austriaca 43.2 24.3 - - - - 
Emys orbicularis 27.8 32.4 11.1 - - 16.5 
Lacerta agilis 31.8 14.4 22.1 - - 21.4 
Lacerta viridis 60.8 19.9 - - - - 
Zootoca vivipara 14.7 = 11.9 = = 57.8 
Natrix natrix 22.9 18.3 19.8 - - 29.0 
Natrix tessellata 48.4 24.2 10.8 - - 10.0 
Vipera berus 16.5 - 17.9 15.0 40.4 - 
Zamenis longissimus 67.6 16.9 = = = = 
Importance for и species: 10 7 7 1 1 6 





optimum suitability is expected to be reached at 28.20 within a narrow range of 28.02 to 
28.32 (using the 10 percentile training presence logistic threshold of 0.4097) (fig. 1). 

Secondly, temperature seasonality (bio 4) was an important contributor to modeling 
potential distribution of 7 species, especially Emys orbicularis (32.4 %). Response curves for 
this variable also displayed a similar bell-shaped response. Both these predictors express 
the homogeneity of temperatures throughout the year, meaning that apparentley most spe- 
cies try to avoid extreme temperature oscillations. 
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Fig. 1. Response curve of isothermality (bio 3), used for predicting the potential distribution of Zamenis longis- 
simus at the local scale (Western Podillya in Ukraine). 
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Fig. 2. Response curve of the Human Footprint (HF), used for predicting the potential distribution of Zootoca 
vivipara at the local scale (Western Podillya in Ukraine). 
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Next to temperature seasonality is the importance of the minimum temperature of the 
coldest month (bio 6). This factor is obviously related to winter survival, and for 7 species 
the percent contribution of the minimum temperature of the coldest month ranges between 
10.8 (for Natrix tessellata) and 22.1 % (for Lacerta agilis). 

Of the non-bioclimatic variables a profound contribution in predicting the species 
geographic distribution models for 6 species has been displayed by the Human Footprint. 
The HF percent contribution reaches values ranging from 10.0 % (for Natrix tessellata) to 
57.8 % (for Zootoca vivipara). Curiously enough, for these species the HF positively affects 
predicted presence probability, as evidenced by the corresponding response curves. In the 
case of Zootoca vivipara, for example, the response to increasing values of the HF resembles 
a sigmoid curve (fig. 2). 

Previous studies have demonstrated that anthropogenic factors have the capacity to 
affect natural species populations (Benayas et al., 2006; French et al., 2010), mostly in a 
negative way. The Human Footprint used in the modeling includes multiple factors related 
to human disturbance (e. g. human population size, land use, infrastructure and the 
degree of human access). One assumption is that these human effects greatly modify and 
generate new environments that are favorable for certain species. Another assumption, 
as evidenced by recent studies, can be linked to net primary productivity (NPP) (Luck, 
2007). It appears that humans preferentially settle in areas of mid to high productivity. 
This pattern is consistent with relationships between NPP and species richness, where it 
appears that richness is typically maximal at high or midrange values of NPP. Hence, it is 
not surprising that positive correlations between human population density and species 
richness for various taxonomic groups have been recorded (for instance, Araújo, 2003). 


Potenial distributions of the species 

Based on the maximum entropy modeling algorithm and using 9 environmental 
variables (VIF < 10), we obtained 10 raster outputs modeling the distribution of the 
considered species. According to table 3, all the species distribution models were better 
than random (AUC > 0.5). For most of the species, models performances were “good” 
(AUC > 0.8) and “excellent” (AUC > 0.9). For only one species, Vipera berus, the model 
performance was “fair” (AUC > 0.7), apparently because of the small number of occurrences 
(п = 6). 

In terms of the percent of suitable habitat in relation to the study area (table 3) the 
species can be ranked in a descending order of nature conservation priority: Zamenis 
longissimus > Lacerta viridis > Zootoca vivipara > Natrix tessellata > Coronella austriaca > 
Emys orbicularis > Anguis fragilis > Natrix natrix > Lacerta agilis. We do not consider 
Vipera berus because of the relatively poor model performance. 


Table 3.Summary statistics for Maxent habitat suitability models 





Species Suitable habitat, %* AUC + SE** 

Anguis fragilis 25.8 0.8840 + 0.006 
Coronella austriaca 22.1 0.9144 + 0.004 
Emys orbicularis 23.0 0.9350 + 0.004 
Lacerta agilis 37.5 0.8623 + 0.004 
Lacerta viridis 9.4 0.9652 + 0.002 
Zootoca vivipara 13.7 0.9426 + 0.005 
Natrix natrix 32.9 0.8898 + 0.004 
Natrix tessellata 15.2 0.9435 + 0.002 
Vipera berus 44.3 0.7670 + 0.012 
Zamenis longissimus 4.6 0.9723 + 0.003 





* % of study area > 10 percentile training presence logistic threshold. 
** SE — standard error. 
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Fig. 3. Areas (downward diagonal filled polygons) in Western Podillya (Ukraine), where the average predicted 
habitat suitability for reptile species exceeds 0.5 (Districts: 1 — Terebovlianskyi, 2 —Husiatynskyi, 3 — Bu- 
chatskyi, 4 — Chortkivskyi, 5 — Chemerovetskyi, 6 — Horodenkivskyi, 7 — Zalishchytskyi, 8 — Borshchivskyi, 
9 — Kamianets-Podilskyi, 10 — Zastavnivskyi, 11 — Khotynskyi). 





Fig. 4. Areas (polygons) in Western Podillya (Ukraine), where there is a predicted probability for the accom- 
modation 9, 8 or 7 reptile species (gradient from dark gray — 9 species to light — 7 species). Districts numbered 
as in fig. 3. 


We combined the predictions into an “ensemble” map by calculating the average of the 
predicted habitat suitability in each pixel for the considered species using a threshold of 0.5. 
Contours enclosing such areas (polygons) are depicted in fig. 3. These areas are likely very 
suitable for the species in general and suggest that surveys should be conducted prior to 
management actions to determine whether any of the species are present and the degree to 
which they may be impacted. In this respect relatively large and promising areas are located 
in Zalishchytskyi, Borshchivskyi and Zastavnetskyi Districts’, especially alongside the 
Dnister, within Chortkivskyi District and around the historic city of Kamianets-Podilskyi. 

By overlaying the binary maps derived for separate species that indicate either presence 
or absence, a summarizing species richness map was produced (fig. 4), emphasizing areas 
(polygons), where there is a predicted probability for accommodating 9, 8 or 7 reptile species. 
As seen, assumed areas of exceptional reptile species richness are located primarily in 
Zalishchytskyi and Zastavnetskyi Districts, the largest one (with the centroid in 25.740064° 
E, 48.648091°N) sharing the territory of both administrative units to the north and south of 
the river Dnister. Centroids of the other eight polygons, assumed to accommodate 9 species, 


2 A rayon (also raion) is a type of administrative unit of several post-Soviet states (such as part of an 
oblast). 
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are: 25.950000° Е, 48.495847° М; 25.845833° Е, 48.508347° М; 25.858333° Е, 48.525014° 
N; 25.829166°Е, 48.541681° М; 25.979166° Е, 48.554181° М; 25.702206° E, 48.763984° М; 
25.733333° E, 48.845847° М. Given the lack of complete information about the distribution 
of the considered species, these locations could be the first to be surveyed for selecting 
potential conservation areas and also be used to verify by ground proof conclusions made 
in the course of the modeling. We also suggest a focus should be made on patches located 
alongside the Zbruch river and the canyon area of Kamianets-Podilskyi. Such informed 
conservation measures should improve the effectiveness of reserve networks and, ultimately, 
will contribute to better protection for biodiversity at the local level. 


Conclusion 


Additional information on the distribution of reptiles (as, for instance, Vipera berus), 
the current status of populations, and the natural history of species in the Western Podillya 
in Ukraine are necessary to develop a conservation program with specific objectives, so 
that guided decisions can help to mitigate negative effects in the populations. Furthermore, 
the protected areas and municipalities in the Western Podillya in Ukraine must develop 
monitoring plans in their areas that contain detailed information on the presence 
or absence of the species. In this respect SDMs can connect species occurrence with 
environmental conditions and be useful tool. SDMs could also become a tool to investigate 
empirical relationships between species occurrences and anthropogenic factors, which may 
be ambiguous or even benefit certain species. Finally, SDMs can be used to identify sets 
of conservation areas, where habitat suitability and factors enhancing species richness are 
in place. Improving data availability and comparing predicted distributions with ground- 
truthing results obtained by undertaking accordant field investigations are necessary for 
future progress in this field of informed biodiversity conservatio 


The authors deeply thank O. Vikyrchak and M. Drebet for their help in collecting the occurrence data. 
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